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Abstract. This paper presents a novel method of approximating the scalar 
Wiener-Hopf equation; and therefore constructing an approximate solution. 
The advantages of this method over the existing methods are reliability and 
explicit error bounds. Additionally the degrees of the polynomials in the ra- 
tional approximation are considerably smaller than in other approaches. 

The need for a numerical solution is motivated by difficulties in computation 
of the exact solution. The approximation developed in this paper is with a 
view of generalisation to matrix Wiener-Hopf for which the exact solution, in 
general, is not known. 

The first part of the paper develops error bounds in Lp for 1 < p < oo. 
These indicate how accurately the solution is approximated in terms of how 
accurate the equation is approximated. 

The second part of the paper describes the approach of approximately 
solving the Wiener-Hopf equation that employs the Rational Caratheodory- 
Fejer Approximation. The method is adapted by constructing a mapping of 
the real line to the unit interval. Numerical examples to demonstrate the use 
of the proposed technique are included (performed on Chebfun), yielding error 
as small as 10~^^ on the whole real line. 



1. Introduction 

The Wiener-Hopf method is used for a broad collection of PDEs which arise in 
acoustic, finance, hydrodynamic, elasticity, potential and electromagnetic theories 
[1, 15]. It is an elegant method based on the exploitation of the analyticity properties 
of the functions. For the scalar Wiener-Hopf the solution can be expressed in terms 
of a Cauchy type integral [16, Ch. 1.3]. 

In more complicated scalar Wiener-Hopf equations the exact solution is difficult or 
slow to compute, see e.g. [2,14,24]. Approximate solutions were considered early 
on but they were mainly constructed using ad hoc observations [16, Ch. 4.5]. In 
2000 a systematic way of approximating the Wiener-Hopf equations was published 
by I. D. Abrahams [3]. Since then it proved popular and found applications in 
different branches of mathematics including finance [10]. 

The method proposed in [3] is based on uniform approximations of the kernel on the 
whole strip by a two point Fade approximation (the two points being and cxi on the 
real axis). However, there are two issues that make the application of this method 
difficult. The first is that it is unclear when the two point Fade approximation 
has small error on the whole strip [6]. Secondly, even if the maximum error of 
approximating the kernel on the whole strip is known, there are no error bounds 
presented to calculate the resulting error in the solution. 
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Another motivation for the development of approximate methods is the matrix 
Wiener-Hopf problem. The determination of a good numerical solution is important 
for the matrix Wiener-Hopf problem, since as yet there is no constructive way of 

solving it in general. 

This paper aims to present a consistent method of approximately solving the scalar 
Wiener-Hopf. The following questions are addressed: 

Problem 1 . Given a scalar Wiener-Hopf problem find an approximate solution which 
can be demonstrated to be within a given accuracy from the exact solution. 

Problem 2. How to perform Problem 1 computationally in a way that is reliable, 
optimal and easy. 

Problem 1 will be addressed in Section 3 with new estimations for the error in 
the Wiener-Hopf factors. This is done though expressing the factors in terms of 
the Hilbert transform, which is connected to the Cauchy singular integral. Note 
that although the method proposed in this paper involves the construction of an 
approximate solution to the given Wiener-Hopf equation, the solution is the exact 
solution of a perturbed Wiener-Hopf equation. 

The main difficulty in Problem 2 is that the Wiener-Hopf equation is set on an 

unbounded interval. We propose the novel method based on an appropriate map- 
ping to the unit interval, see Section 4. This allows the Caratheodory-Fejer rational 
approximation to then be used in Section 5. Lastly, numerical examples are given. 

2. Preliminaries 

The following conventions will be used throughout the paper. A strip around the 
real axis (given t_ < < t+) is {a = a + iT : T- < t < r+}c C. The subscript -|- 
(or — ) indicates that the function is analytic in the half-plane t > t- (or r < r+). 
Functions in the Wiener-Hopf equation without the subscript are analytic in the 
strip. 

The Wiener-Hopf problems are recalled below. 

The multiplicative Wiener-Hopf problem is: given a function K (analytic, zero-free 
and (a) —> 1 as I cr I oo in the strip) to find functions and K- which satisfy 
the following equation in the strip: 

K{a) = K_{a)K+{a). (2.1) 

In addition K- and are required to be: 

• analytic and non-zero in the respective half-plane. 

• of subexponential growth the respective half-planes: 

I logii'±(a)| = 0(|a|^), p < 1, as |<t| ^- oo. 

The function K{a) in the multiplicative Wiener-Hopf factorisation will be referred 
to as the kernel. The functions and will be called factors of K. Those factor 
are unique up to a constant. In other words if there are two such factorisations 
K = K+K_ and K = P+P_ then: 

K+ = cP+ and K_=c~^P-, 

where c is an arbitrary complex number [13]. In this paper they will be normalised 
so that K±{a) ^ 1 as |cr| — >■ oo in the strip. 
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Remark 2.1. The above normalisation at infinity allows to control the regularity of 
the Fourier transform of an approximation. This is useful if a factorisation of the 
Fourier transform is used, say, for a solution of differential equations. 

The multiplicative Wiener-Hopf factorisation can be reduced to the additive Wiener- 
Hopf problem via application of logarithm. Under the conditions on the kernel K, 
let / = log if then additive Wiener-Hopf problem is: given a function / to find 
functions /+ and /_ which satisfy the following equation in the strip: 

/(a) = /_(a) + /+(a). (2.2) 

Additive and multiplicative problems are equivalent only in the scalar case. In the 
matrix Wiener-Hopf and K- will in general not be commutative so no such 
equivalence is possible. 

The existence of solutions to the additive and multiplicative scalar Wiener-Hopf 
problem is addressed in for example [16, Ch. 1.3]. 

These two Wiener-Hopf splittings are key to solving the general Wiener-Hopf prob- 
lem. That is: given functions A and C find functions and which satisfy 
the following equation in the strip: 

A{a)^+ (a) -I- (a) -|- C{a) =0. (2.3) 

For more details about the Wiener-Hopf problem see [16] . 

2.1. The Riemann-Hilbert Problem. The Wiener-Hopf problem is a special 
case of the Riemann-Hilbert problem. Roughly, the Riemann-Hilbert problem has 

a more general contour instead of the strip and the conditions on the function A{a) 
are weakened. In particular it is required that A[a) is only Holder continuous on 
the contour. 

Let G{t), g{t) be Holder continuous functions on a simple contour S. The Riemann- 
Hilbert problem is to find an analytic function E{a) which has values E^{t), E^{t) 
on the contour as the limit is taken from different sides of the contour and which 
satisfies the equation: 

G{t)E_{t) + E+{t) -\- g{t) = 0, t e S. 

The Wiener-Hopf equation (2.3) can be considered as the special case of the Riemann- 
Hilbert problem so a good approximation of the values of the function on the real 
line is sufficient. This is simpler than trying to construct an approximation that 
agrees well on the whole strip of analyticity. 

3. Estimates on the Approximate Wiener-Hopf Factorisation 

The approximate solution to (2.1) can be foimd by approximating the kernel in the 
Wiener-Hopf equation. This section quantifies the difference between the modified 
and the original Wiener-Hopf equation, to address Problem 1 introduced earlier. 
More precisely, let \K{a) — K{a)\p < ep, the aim is to bound \K±{a) — K±{a)\p in 
terms of Cp and computable quantities of K{a). 

The first step is to link the Wiener-Hopf factorisation to the Hilbert Transform 
H{f){y) of/withyeM. If: 

m = f+{y) + f-{y), 

then [17, Ch. 2]: 

iH{f){y) = U{y) - f_{y), 
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and hence: 

f±{y) = lfiy)±'-H{f){y). 

The advantage of expressing f± in terms of the Hilbert transform is the following 
theorem. 

Theorem 3.1 (Titchmarsh-Riesz). [17, p. 92] Let f e Lp{R) for some 1 <p <oo, 
then: 

\\H{f)\\,<c{p)\\f\\„ 

where the best constant is: 

/ tan(7r/(2p)) if l<p<2, 
""^'-l cot(7r/(2p)) if 2<p<oo. ^'^ > 

Based on these classical estimations we obtain the following bounds on the additive 
Wiener-Hopf factorisation. 

Lemma 3.2 (Additive Bounds in Lp ioi 1 < p < oo). Let f{y) = f+{y) + f-{y) 
and f{y) = f+{y) + with \ \f{y) - f{y)\\p < Cp then: 

\\f±{y)-f±{y)\\p<l{i + c{p))ep, 

where c{p) is defined as in (3.1). 

Proof. Expressing in terms of the Hilbert transform and using Minkowski's inequal- 
ity: 

||/+(y)-/+(y)||p = \\^f{y)+'-Hif){y)-^f{y)-'-H{f){y)\\p 

< wlm - Imwp + \\'^H{f)iy) - ^i/(/)(2/)iip 

< \{^ + c{p))ep. 

Here the linearity property of the Hilbert transform, H{f) + H{g) = H{f + g) is 
used. □ 

The theorem gives bounds on the error on the real line but in fact since the functions 
are analytic, the bounds will hold in the upper/lower half-planes by the maximum 
modulus principle. 

It will be useful to express the Wiener-Hopf factors in terms of the Hilbert trans- 
form; 

K±{y) = exp {\\ogK{y) ± '-H{\ogK){y)) = K^/\y)e^p {± '-HilogK){y)). 
The main new result of this section is: 

Theorem 3.3 (Multiplicative Bounds in Lp for 1 < p < oo). Let K{y) and K{y) 
be two kernels and m < \\K\\p < M. If \\K{y) — K{y)\\p < Cp then: 

\\K^{y) - KM\p < ^^^ 2(m - O + '^P^^'-' 



where c{p) is defined in (3.1). 
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Proof. The first step is to express ttie product factorisation K{y) = K+{y)K-{y) 
as an additive factorisation. Taking logarithms: 

log K{y) = logK+{y) + logK_{y). 

From the assumptions on the kernel K{y), \ogK{y) will be a continuous and single 
valued function and hence the additive decomposition is well defined. Given: 

\\K{y)-k{y)\\^<e^. 
Then using the mean value inequality (and ||.^(?/)||p > m — Cp): 

II \og{K{y)) - \og{k{y))\\, < ^^e^. 

Apply Lemma 3.2, which gives: 

\\\og{K^{y)) - log(^±(2/))||p < J (1 + c(p))ep. 

Lastly, note that \\k±{y)\\p < (M + ej,)^^ exp(£^) since: 



\\K±{y)\\, = \\K'/\y)e^p{±'-Hilogk){y))\\, 

< (M + ep)i/2|| cxp ( ± '-Hi\ogK){y))\\p 

< (M + ep)V2|| exp ( ± ^<^{H {log k){y))% 



< (M + e,)V2|| exp ( ± i(//(Arg(^))(y)))||, 

< (M + e,)V2exp(^), 



and so: 

\\K±{y) - k±{y)\\, < (M + ep)i/2exp(^)||log(if±(2/))-log(if±(t/))||p, 

from which the result follows. □ 

Remark 3.4 (Real Kernels). In the case when kernels K and K are real valued, the 
bounds could be simplified, since then ||.ftr±(t/)||p < (M + e^)^/^ and so: 

(M-\-f ^1/2 

\\K±{y) - K±{y)\\, < i— pA-(l + c{p))e,. 

Z\Tfl tp) 

Remark 3.5 (Loo norm). The above theorem does not include the case of L^o norm, 
in fact no such result is possible for the Loo 

norm as is demonstrated by the counter 

example below. 

Consider Ki{y) = 1 defined on the real line and 

{1 + e ifO<j/< arccot(ni); 

1 + e - if arccot(ni) <y< arccot(n2); 

1 if arccot(n2) < y < oo, 

and 

K2{y) = K2{-y). 
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Then: 

\\K2{v)-KMoo=e, 
but it has been shown by making n2 — ni small the factors can differ by an arbitrary 
large amount in Lqo norm, see [5]. 

Below, the proposed method of the solution of Problem 1 is presented: 

• Approximate, with arbitrary accuracy, the Wiener-Hopf kernel K by ra- 
tional functions K, that can be done by Runge's approximation theorem [18] . 

• Perform the Wiener-Hopf factorisation of the rational kernel K by inspec- 
tion. The K+ will have the form: 

(y-zi)...(y-z„) 

(y-pi) •■•(?;- Pm)' 

where the Zi (pj) are all the zeroes (poles) of K which lie in the lower half- 
plane. The other factor K- will have the remainder poles and zeroes in the 
upper half-plane. 

• Using Theorem 3.3 the error — K±{y)\\p can be calculated. 

The rest of the paper will concentrate on Problem 2 i.e. the practical aspect of 
Problem 1. 

4. Mapping of the Real Line 

In the previous section it was proved that the approximation of the kernel results in 
a computable error in the Wiener-Hopf factors. Thus, to simplify the problem one 
may wish to approximate a kernel by another which is easier to factorise e.g. by 
rational functions. To construct such an approximation the first step is to employ 
a mapping of the real line to the unit circle or the interval [—1,1]. Those mappings 
transform the problem of approximation on the whole real line to a simpler problem 
on the unit circle or the interval. Methods for approximating on the unit circle and 
the interval will be discussed in the next section. 

At this point it is important to fix notation. There will be functions defined on 
the real axis, the unit circle, and the interval [—1, 1]. To easily distinguish between 
them the following conventions will be used: 

• Functions on the real line will be denoted by capital letter and variable y, 

for example F(y); 

• Functions on the unit circle will be denoted by small case letter and variable 

z, for example f{z); 

• Functions on the interval [—1,1] will be denoted by bold capital letters and 

variable x, for example F(x); 

• The rational approximation for functions will be denoted by letters with a 
tilde, for example F{y). 

It is typical in applications for the kernel to be an even function (i.e. F{y) = F{—y)) 
and the mapping is simplified in this case. Therefore it is instructive to consider 
the mapping for even functions before considering the general case. 
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4.1. Even Functions. Define a conformal Mobious map M{y) that maps the real 
line to the unit circle as follows [7, Ch 13]: 

1 i, — 1 1-)- —i, oo I-)- 1. 

The easiest way to picture the map is to think of the plane plus infinity as a 
Riemann Sphere and the map as a composition of two rotations of the Riemann 
Sphere. Explicitly it is given by: 

The next map is a projection of the unit circle to the interval [—1,1] by the 
Joukowski map (a conformal map which, incidentally, comes from aerodynamics 
[20]). The Joukowski map is: 

J{z) = liz + z-'). (4.1) 

Restricted to the unit circle the map returns the real part of z. This step requires 

the function to be even. Even functions will be mapped to a function on the unit 
circle with the property f{z) = f{z~-^). Composing these maps together: 



This enables the construction of a function F(a;) on the interval from a given func- 
tion F{y) on the real line as follows: 

F(.t) = FiM-^J~\x)). 

Then the algorithm of the next section can be used to rationally approximate F{x) 
on [—1, 1] by F(a;). Suppose that the error in this approximation in the Lqo is A. 

Then map the function back to the real line by: 

F{y) = F{JM{y)). 

Furthermore F{y) is a rational function on the real line that approximates F{y) on 
the real line with the L^o error at most A. 

The maps described above are combined with the Carathcodory-Fejcr algorithm in 
Chebfun to give the code used for the numerical examples in the next section. 

Remark 4.1. If the kernel is odd i.e. F{y) = —F{—y) the kernel can be squared 
to produce an even kernel. Then the approximation can be produced as above and 
square rooted. Note that the explicit multiplicative Wiener-Hopf factorisation of 
the square root of a rational function is done by inspection in the same way as for 
the rational function. Also note that any function can be written as a sum of even 
and odd function, thus the above extends to a general function. 

4.2. General Functions. For a general kernel the above map can be modified. 

The above method contains the following ideas. The real line is conformally mapped 
to the circle. Then only half of the circle is mapped to the interval [—1, 1]; however 
this poses no problem since the other half is the same due to the function being 
even. The function is approximated and mapped back to the unit circle. The 
modification is needed to ensure that half of the function is not discarded. Thus, 
after mapping the real line to the circle, apply the map S : z ^ z^; this makes the 
values run twice as fast and the whole function now fits on the half circle. Next, 
as before, map this half circle to the interval, approximate, and map back. The 
values that are taken on the half circle are spread out to the circle by the map 
: z ^ \fz. Calculating the resulting map gives: 
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Figure 1. Showing the map JSM(y) of the real Une (plotted from 
-20 to 20) to the interval [-1, 1]. 



JSM{y) = ^^, M-'S-'j-\x)^^=^, (4.2) 

see Figure 1. 

Remark 4.2. Before using a map it is beneficial to apply a rescaling and a shift 
(i.e. Mobious map) to the real line so that the part of the function which has a 
fast changing gradient fits in the interval [—5, 5]. This part of the real line is well 
resolved, see Figure 1. 

4.3. Associated Orthogonal Rational Functions. In this section the proper- 
ties of the map (4.2) will be studied. This will be done by considering the new basis 
functions which are created by this map. A starting point is to choose a basis func- 
tion for the circle; a good choice for this is the Fourier basis. Then the Joukowski 
map will create the new basis for the interval [—1,1] as the image of the Fourier 
basis, under the change of coordinates, [9, Sec 1.6]. This new basis is the Chebyshev 
polynomials of the first kind (they are also the Fabe polynomials on [—1,1] up to a 
multiplicative factor [19]). Furthermore, the mapping M"^ S^^ J^^ {x) from [—1, 1] 
to [— c»,c»] will create a new basis (from Chebyshev polynomials) called Rational 
Chebyshev Functions, TBn [8]. They are defined as: 

TS„(y) = T„(- y ^ 



where the T„ are the Chebyshev polynomials of the first kind and . ^ is the 



uals of the 

map JSM{y), see (4.2). Note that, despite the name, TBn are not all rational: 
the TB2n+i are rational functions divided by ^ (y^ -|- l) (the TB2n are all rational 
function) . An important property of TBn is orthogonality: 

TBn{y)TB„,{y) ^ / 7r/2 for n = m > 0; 



J 1 + 1 ioT n ^ m. 

This is derived from the corresponding orthogonality property of the Chebyshev 
polynomials. This property can be used to calculate the coefficients in the expansion 
of a function in terms of the Rational Chebyshev Functions. 
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Another useful property of TBn is the shape of the domain of convergence. This 
is especially informative when the kernel has no singularities on the whole real line 
including infinity. For standard basis functions the domain of convergence has a 
specific shape and the size is determined by the position of singularities [8]. For 
example, Taylor series converge on a circle and Chebyshev series converge on an 
ellipse, the sizes of which are dictated by the position of the singularities. The 
shape of the domain for TB^ is the exterior of two circles which are the image of 
lines parallel to the real axis under the Mobious map '^^j^J . This is derived from 
the corresponding domain of convergence for Chebyshev polynomials which in turn 
comes from the Fourier scries [8]. 

So in particular those base functions are well suited for functions which have singu- 
larities in the bounded part of the complex plane, like Example 6.1. This compares 
favourably to the domain of convergence of Taylor expansion for a function like the 
one in Example 6.1. In particular, this explains why all methods of approximation 
which are based on the Taylor expansion do not produce good results on the whole 
real line (e.g. one point Fade approximation). 

For functions with a singularity at infinity the following theorem by Boyd is applic- 
able. 

Theorem 4.3 (Convergence of the Fourier series for an algebraically decaying 
or asymptotically constant function). [8, Sec 5] // a function u{y) is free from 
singularities on the real axis and has the inverse power expansion 

w ~ Co + ci/y + C2/y^ + ... (4.3) 

as y ^ oo and a similar series asy ^ —oo, then the coefficients of its representation 
as a Fourier series in the new coordinate t (where z = e^'^'* and z = SM{y), see 

(4.2);, 

oo oo 

u{y) = ^ ctn cos(ni) + '^b„ sin(nt), (4.4) 

n—Q n—1 

will have exponential convergence in the sense that |a„| and decrease with n 
faster than any finite inverse power of n. 

In the case of (4.3) converging then (4.4) converges geometrically. 

In the case of (4.3) being asymptotic hut divergent (then u{y) will he singular at 
infinity hut all the derivatives will be hounded there) then (4.4) converges suhgeo- 
metrically. 

Remark 4.4. The theorem covers precisely the class of functions that arises as 
kernels in the additive and multiplicative Wiener-Hopf factorisation. In the additive 
case the kernel must decrease faster than some power and in the multiplicative case 
they are asymptotically 1. 

The above theorem shows one of the strengths of the TBn compared to other base 
functions on an unbounded interval e.g. the Hermite and sine functions. The latter 
only converge algebraically for algebraically decaying functions. 

5. Rational Approximation 

The best rational approximation in Loo is not practical to use, see [23] . This section 
describes a near best method of rational approximation named AAK (Adamjan- 
Arov-Krem) [4] in the general case and Rational Caratheodory-Fejer [11, 21] in 
the case of real valued functions. The AAK theory is briefly reviewed at the end. 
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Unfortunately it seems that there is no software that performs AAK approximation; 
only software for Rational Caratheodory-Fejer have been developed that is Chebfun 
(MATLAB) [23]. 

5.1. Real Rational Approximation. Below a brief summary of the method from 
[22] is presented. Let F{x) be a real continuous function on [—1, 1]. For any natural 
number M, F{x) has a partial Chebyshev expansion: 

M 

F{x) = ^akTk{x) + Gm{x) = Fm{x) + Gm{x), 

where 

2 c\r 
ak = - F{x)Tk{x) 
^ J-i 



and Tfe (x) arc Chebyshev polynomials of the first kind. 

An application of the Joukowski map (see (4.1)) to Fm{x) produces an associated 
function on the unit circle: 

1 ^ 
Fm{x) = -fM{z), where Jm{z)= ^ amZ^ . 

k=-M 

The above is true because of the intimate connection of the Chebyshev polynomials 
and Joukowski transform given by: 

Tk{x) = ^iz'' + z-^). 

Define: 

M 



k—m—n 

where we are seeking a [m,n] rational approximation. Next approximate f~^{z) on 
the unit circle by an extended rational function of the form: 

r{z) = - \ , (5.1) 

2^k=o ^kz 

where the numerator is a bounded analytic function in [zl > 1 and the denominator 
has no zeroes in |2;| < 1. Call a class of functions of the above form Rm,n- Note 
that this is not a typical class to approximate by, but it is the one for which a neat 
solution exists. After the approximation (5.1) is obtained it can be truncated to 
get a rational approximation. To obtain an approximation in Rm,n consider a real 
symmetric Hankel matrix: 

^ ^m— n+l • • • \ 

x= . 

\ am ... J 

Let Xi be eigenvalues of X arranged in decreasing order by magnitude of absolute 
value. And let (ui, . . . ,UM+n-m) be an eigenvector for the A„+i eigenvalue. Then 
the following theorem is true: 

Theorem 5.1 (Takagi). The analytic function /+ has a unique best approximation 
f on the unit circle = 1 in Rm,n given by: 

f{z) = b{z) - f+{z) 

where b is the finite Blaschke product 
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M 'J^i "I" ■ ■ ■ "I" UM+n—m' 
UM+n-m + • • • + UiZ 



.M+n—m 



b{z) = Xn+lZ 



■ M+n—m—l ' 



And the approximation error in L^o norm on the unit circle is |A„+i|. 

Once the approximation is obtained, it is mapped back from the unit circle to the 
unit interval by the inverse of the Joukowski transform. A subsequent truncation 
gives a near best rational approximation on the interval. 

This is a fast and efficient way of computing rational approximations. Furthermore, 
it enables to predict the error on the approximation even before the approximation 
is computed by considering eigenvalues. 

5.2. AAK Approximation. This section presents more general results to the 
ones covered in the previous section. These apply to complex valued functions on 
the unit circle and are taken from [4]. 

Definition 1 . Given a natural number n, define a class of bounded functions 
on the unit circle which can be expressed as: 



where r{z) is a rational function that has no more than n poles all inside the unit 
circle and h{z) € Hoo- 

The AAK approximation solves the following problem: 

Given a function f{z) € and n natural number, find the best approximation in 
the Lao norm from functions in the H^. 

In other words it is required to find: 



The solution to the above is presented in the next theorem. 

Theorem 5.2. [4] Let f{z) € L^o, then D„{f) = .s„(M) where M is the Hankel 
matrix build out of Fourier coefficients of f{z) and s„ is the n*'' singular value of 
it. Moreover, 



g{z) = r{z) + h{z), 





OO 



i=i 



The next section illustrates theory with examples. 



6. Numerical Examples 



In this section the numerical examples are provided to illustrate the theory. All the 
numerical experiments are performed in Chebfun (MATLAB). 
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Figure 2. The function to be approximated F{y) — (^5^- 
6.1. Example la. The first example (see Figure 2) is 



Take the finite branch cut from i to ki and from —i to —ki. This kernel is closely 
associated with the matrix kernel factorisation from problems in acoustics and 
elasticity and was studied by I. D Abrahams [3]. The case k = 2 will be considered 
in the numerical examples. The first approach is to approximate this given function 
on an interval, a process called domain truncation. On the surface, the results 
look promising (see Figure 3); the error decreases as the degree of the numerator 
polynomial is increased. Nevertheless there are problems, the most obvious being 
that [20,4] or any other approximations plotted in Figure 3 have very different 
behaviours at infinity than F(jj). The error will be small but only on the given 
interval and will not be controlled outside it. It might be tempting to try to rectify 
this, since F{y) — 1 will be zero outside a sufficiently large interval. One might 
suggest to look at [n,m] approximations where n is smaller than m: this type of 
function will go to zero at infinity. But, although the error on most of the real axis 
is very small, just outside the interval it tends to be greater (and much larger than 
it is inside the interval) . Another feature of this type of approximation is how the 
poles and zeroes are positioned. Although some of the poles are positioned on the 
branch cuts [i, 2i] and [—i,—2i], there are a lot of singularities introduced elsewhere 
(see Figure 4). 

In the next section, the mapping constructed in Section 4 will be used to overcome 
the described problems. This will also result in smaller degree of the polynomials. 

6.2. Example lb. The previous example (6.1) is treated with the mapping pro- 
posed in Section 4 for even functions. The function mapped to the interval [—1, 1] 
becomes: 

In Figure 5 the new function is plotted; notice that in this form it is much easier 
to approximate, since it has a slow changing gradient. This is confirmed by the 
error curve for the [4,4] approximation; the error is of order 10^^. Furthermore, 
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Figure 3. The error of [20,4] (blue) [40,4] (red) [60,4] (green) 
rational approximations calculated using Chebfun k = 2. 
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Figure 4. The poles and zeroes of [20,4] [40,4] [60,4] rational approximations. 



a better result is obtained with [5, 5] approximation with error of order 10^^^ (see 
Figure 6). At this point the approximated function can be mapped back to the 
real line and the Loo norm of the error will stay the same, see Figure 7 (solid green 
curves). Note that it now becomes a [10,10] rational approximation. This example 
was chosen because the Wiener- Hopf factors can be easily seen by inspection. They 




Figure 6. The error in approximating F(a;) by [5,5 



are: 

The factors of the rational approximation can be easily calculated and compared to 
the exact solution above; this is plotted in Figure 7 (dotted red curve). This error 
is smaller than the error obtained in [3] even with [30, 30] for the same function. 
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Figure 7. The error in approximating F{y) by [10,10] (solid green 
curves) and the error after it is spht in factors, \F+{y) — F+(y)| 
(dotted red curve). 
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Figure 8. The poles and zeroes in approximating F{y) by [10,10] 
the poles are filled in red circles whereas the zeroes are green circles. 
The scale of the x-axis is 10""^. 
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What is more, the poles and zeroes also lie almost exactly on the branch cuts [i, 2i] 
and [—i,—2i\. The close proximity of zeroes and poles mimic the behaviour of a 
branch cut. This indicates that the behaviour in the whole complex plane is correct 

(sec Figure 8). Note that the scale of the x-axis is 10"''. 

Bounds developed for real-valued kernels in Section 3 could be applied to this 
example. The L2 bounds will be used. Note c(2) = 1, M = 1 and m = 1/2 this 
gives: 

\\Fi{y)-F±{y)\\2<2\\F{y)-F{y)\\2. 

The 1/2 norm is easy to estimate numerically in Chebfun by using the norm com- 
mand. For the [10,10] rational approximation these calculations give: 

-%)||2 = 9.5 X 10-10, and \\F±{y) - F±{y)\\2 = 6.5 x IQ-^'', 
which agrees with the bounds above. 

6.3. Example 2. The function is (Figure 9): 



i^(y) = ^^^tanh,. 

This kernel is related to the one considered in Koiter's paper of 1954 [12] and is 

typical for electrostatic and slow flow problems. In the original article, rational 
approximations were considered by choosing the coefficients by hand, accuracy of 
10"^ was achieved. Once again the function becomes easier to approximate (on the 
whole real line) once it is mapped to the interval [—1,1], Figure 9. The mapped 
function is (the general mapping is used): 

1 'T- 
F{x) = — tanh 



X VI - a;2 ■ 

This can then be approximated with errors of 10^^ (giving the graph similar to 
Figure 6) with [12, 12] even though the point at infinity is singular (it is an accu- 
mulation point of a sequence of poles). Interestingly, the position of zeroes and 
poles (see Figure 10) is very similar to that found in [3, Fig 8] with a different 
method of approximating. This kernel has exact factorisation which is given by: 

V^r(l-zy/7r) ^ ' 

But MATLAB does not have an inbuilt complex F function hence the results of 
this approximation cannot be compared to the exact factorisation. Even if there 
was an inbuilt complex F function to compare it to the approximate solution the 
inbuilt function would need to have very high accuracy. This demonstrate that if 
the numerical values of the solution is needed the approximate solution is as good 
as the exact solution. 



7. Conclusion 

This paper demonstrates with use of examples the method of approximate solution 
of the Wiener-Hopf equation. The bounds presented allow to predict the error in 
the approximation. The degrees of polynomials in the rational approximation are 
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Figure 9. The initial function to be approximated F{y) = 
tanhy (top) and the mapped function to be approximated 



F(a;) = i tanh 



(bottom). 



smah which allow to simplify the initial problem significantly. The next project 
would be to generalise the approach to the matrix Wiener-Hopf problem. 
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